% this function is use for generation the params Cn2(z)
% which is based on the HVB model; refer in p23 <<principles of adaptive
% optics>> 3ed
% Written by Lu benchu, Taiyuan University of Technology
% email: benchuul@163.com
function Cn2 = HVB_Cn2(params_profile)
% OUTPUT:
% Cn2(z) function handle with param 'z'which is unit in 'km'.

% Checking the default params setting
if (isfield(params_profile,'A') && isfield(params_profile,'W'))
    A = params_profile.A;
    W = params_profile.W;
elseif (isfield(params_profile,'lambda') && ...
        isfield(params_profile,'theta_0') && isfield(params_profile,'r0'))
    lambda = params_profile.lambda * 10^(6); % [m] to [μm]
    theta_0 = params_profile.theta_0 * 10^(6); % [rad] to [μrad]
    r0 = params_profile.r0 * 10^(2); % [m] to [cm]
    A = HVB_A();
    W = HVB_W();
else
    msg = 'lacking of params in params_profile: A and W or r0, theta_0 and lambda';
    error(msg);
end

if isfield(params_profile,'h_site')
    h_site = params_profile.h_site;
else
    h_site = 0;
end

Cn2 = @(z) 5.94*10^(-23)*z.^(10).*exp(-z)*(W/27)^(2) + ...
    2.7*10^(-16)*exp(-2*z/3) + A*exp(-10*(z-h_site));

    % this sub-function is use for computing the params 'W'
    % ----------------------------------------------------
    function W = HVB_W() % error
        W = 27*(75*theta_0^(-5/3)*lambda^(2) - 0.14)^(1/2);
    % end of the function 'HVB_W'
    end

    % this sub-function is use for computing the params 'A'
    % ----------------------------------------------------
    function A = HVB_A()
        A = 1.29*10^(-12)*r0^(-5/3)*lambda^(2) -...
            1.61*10^(-13)*theta_0^(-5/3)*lambda^(2) -...
            3.89*10^(-15);
    % end of the function 'HVB_A'
    end

% end of the function
end


